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Nmcap-2k  is  a  three-dimensional  computer  code  that  models  interactions  between 
spacecraft  and  plasma  environments  in  low-Earth,  geosynchronous,  auroral,  and 
interplanetary  orbits.  The  code  builds  on  physical  principles,  mathematical  algorithms,  and 
user  experience  developed  over  three  decades  of  spacecraft  charging  research.  Nascap-2k 
has  improved  numeric  techniques,  a  modern  user  interface,  and  a  simple,  interactive  satellite 
surface  definition  module  (Object  Too  IK  it). 

Capabilities  include  surface  charging  in  geosynchronous  and  interplanetary  orbits, 
sbeath  and  wake  structure  and  current  collection  in  low-Earth  orbits,  and  auroral  charging. 
External  potential  structure  and  particle  trajectories  are  computed  using  a  finite  element 
method  on  a  nested  grid  structure  and  may  be  visualized  within  the  Nascap~2k  interface. 
Space  charge  can  be  treated  either  analytically,  self-consistently  with  particle  trajectories,  or 
consistent  with  imported  plume  densities.  Particle-in-cell  (PIC)  capabilities  are  available  to 
study  dynamic  plasma  effects. 

Previously,  we  reported  on  the  accuracy  of  Nascap-2kfs  charging  and  current  collection 
calculations  by  comparing  computed  currents  and  potentials  with  analytic  results,  and  by 
comparing  Nascap~2k  results  with  published  calculations  using  the  earlier  lower  resolution 
codes,  NASCAP/GEO,  NASCA P/LEO,  and  POLAR.  Here  we  examine  the  accuracy  and 
limitations  of  two  new  capabilities  of  Nascap-2k:  modeling  of  plasma  plumes  such  as 
generated  by  electric  thrusters  and  enhanced  PfC  computational  capabilities. 

Nos  cap- 2k  models  one  or  more  ion  engine  plumes  in  full  three-dimensional  geometry, 
including  plume-plume  and  plume-spacecraft  interactions.  The  primary  thruster  beam, 
parameters  describing  the  neutral  efflux,  and  the  initial  charge-exchange  plume  are  imported 
from  a  PlumeTooFgenemted  file.  Nascap-2k  generates  and  tracks  charge-exchange  ions  to  obtain 
plasma  densities  and  calculates  potentials  consistent  with  plasma  densities  and  object  surfaces. 

Nascap-2k's  PIC  capability  has  been  expanded  to  include  boundary  injection,  particle 
splitting,  and  substep  charge  deposition.  The  boundary  injection  replaces  collected  particles 
in  long  running  calculations.  The  particle  splitting  allows  for  modeling  the  effects  of  the 
thermal  distribution  of  velocities,  as  well  as  accommodating  particle  weight  to  variable  grid 
cell  volume.  The  substep  charge  deposition  makes  possible  calculations  in  which  two  effects 
have  significantly  different  timescales.  We  use  calculations  for  simple  geometries  to  explore 
the  accuracy  and  limitations  of  these  capabilities. 
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I.  Introduction 

esigners  of  spacecraft  for  government,  commercial,  and  research  purposes  require  advanced  modeling 
capabilities  to  guide  the  design  of  satellites  that  can  survive  and  operate  property  in  the  natural  and  induced 
environment.  Computer  modeling  of  flight  experiments  (including  SCATHA,  the  SPEAR1 :  series,  and  CHAWS3) 
has  reliably  predicted  both  steady-state  and  dynamic  interactions  between  high-voltage  spacecraft  and  the  ambient 
plasma.  Computer  modeling  has  also  been  applied  to  inherently  dynamic  problems,  such  as  those  involving  three- 
dimensional  space  charge  sheath  formation,  current  flow  in  the  quasi-neutral  presheath,  breakdown  phenomena, 
plasma  kinetics,  ionization  processes,  and  the  effect  of  unsteady  processes  on  spacecraft  charging. 

Nascap-2k  is  being  developed  by  Science  Applications  International  Corporation  (SAIC)  to  address  this  need. 
The  development  has  been  sponsored  jointly  by  the  Air  Force  Research  Laboratory  at  Hanscom  AFB  and  by 
NASA's  Space  Environments  and  Effects  (SEE)  Program  at  Marshall  Space  Flight  Center.  The  current  release  is 
Version  3.  t  and  Version  3,2  is  expected  to  be  released  before  the  end  of  2007, 

Nascap-2k  is  a  widely  used  interactive  toolkit  for  studying  plasma  interactions  with  realistic  spacecraft  in  three 
dimensions.  It  can  model  interactions  that  occur  in  tenuous  (e.g.,  GEO  orbit  or  interplanetary  missions)  and  in  dense 
LEO  orbit)  plasma  environments.  It  incorporates  the  physics  and  numeric  techniques  developed  over  the  last 
thirty  years  to  address  these  interactions. 

Nascap-2k  has  primarily  been  used  to  model  spacecraft  surface  charging  in  tenuous  plasmas,  such  as  at 
geosynchronous  orbit,  and  current  collection  from  dense  plasmas,  such  as  in  low-earth  orbit.  Validation  of  these 
capabilities  has  been  addressed  in  a  previous  publication,4  The  code  can  also  be  used  to  examine  particle  trajectories, 
wake  structure,  and  charging  in  dense  plasmas.  Recently  two  capabilities  have  been  significantly  enhanced.  The 
inclusion  of  ion  densities  from  plasma  sources  (such  as  thruster  plumes)  in  the  computation  of  potentials  in  space 
has  been  enhanced  by  the  direct  input  of  the  ion  plume  and  the  self-consistent  calculation  of  charge  exchange.  The 
particle- in-cel  I  (PIC)  capabilities  have  been  enhanced  to  make  PIC  calculations  more  flexible,  robust,  and  faster 


Figure  I.  Views  of  Nascap-2k  user  interface. 
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//.  Nascap-2k 

Nascap-2k?*  7'  *  *  is  a  three-dimensional  spacecraft  plasma  environment  interactions  computer  code  that 
simluates  a  wide  variety  of  plasma  phenomena.  These  include  spacecraft  charging  in  geosynchronous, 
interplanetary,  auroral,  and  low-Earth-orbit  plasmas,  volume  potentials,  particle  trajectories,  and  resulting  variations 
in  plasma  density. 

Nascap-2k  is  targeted  to  spacecraft  design  engineers,  spacecraft  charging  researchers,  and  aerospace  engineering 
students.  The  graphical  user  interface  is  designed  to  help  less  experienced  users  easily  solve  moderately  complex 
plasma  interactions  problems  white  allowing  the  expert  user  to  tackle  questions  that  have  not  been  previously 
contemplated. 

The  core  capabilities  of  Nascap-2k  are  the  following: 

1.  Define  spacecraft  surfaces  and  geometry  and  the  structure  of  the  computational  space  surrounding  the  spacecraft; 

2.  Solve  for  time-dependent  potentials  on  spacecraft  surfaces; 

3.  Solve  the  electrostatic  potential  about  the  object,  with  flexible  boundary  conditions  on  the  object  and  with  space- 

charge  computed  either  fully  by  particles,  fully  analytically,  or  in  a  hybrid  manner; 

4.  Generate  and  track  electrons  and  ions,  represented  as  macroparti des,  including  computing  the  resulting  surface 

and  volume,  current  and  charge  densities; 

5.  View  surface  potentials,  space  potentials,  particle  trajectories,  and  time-dependent  potentials  and  currents. 

The  earliest  and  most  common  application  of  Nascap~2k  is  to  study  charging  of  spacecraft  in  geostationary  orbit, 
Nascap-2k  is  designed  to  make  this  type  of  analysis  particularly  easy.  Nascap-2k  calculates  surface  potentials  and 
electric  fields  using  the  Boundary  Element  Method  (BEM),10  thus  a  simple  charging  calculation  does  not  require  a 
spatial  grid.  Nascap~2k  contains  a  selection  of  predefined  Maxwellian,  Double  Maxwellian,  and  auroral  plasma 
environments,  which  arc  readily  modified  to  create  custom  environments.  Other  aspects  of  the  environment  include 
the  magnetic  field,  the  sun  direction,  and  the  sun  intensity. 

To  solve  for  the  electrostatic  potential  about  the  object,  Nascap-2k  uses  a  high-order  finite  element 
representation  of  the  electrostatic  potential  that  assures  that  electric  fields  are  strictly  continuous  throughout  space. 
The  electrostatic  potential  solver11  uses  a  conjugate  gradient  technique  to  solve  for  the  potentials  and  fields  on  the 
spacecraft  surface  and  through  the  surrounding  space.  Several  analytic  and  numerical  space  charge  density  models 
are  available,  including  Laplacian,  Linear,  Non-linear,  Frozen  Ions,  Full  Trajectory  Ions,  Hybrid  PIC  (appropriate  to 
the  several  microsecond  timescale  response  to  a  negative  pulse),  and  Full  PIC, 

Particle  tracking  is  used  to  study  sheath  currents,  to  study  detector  response,  to  generate  steady-state  charge 
densities,  or  to  generate  space  charge  evolution  for  dynamic  calculations.  Nascap-2k  generates  macro  particles  (each 
of  which  represents  a  collection  of  particles)  either  at  a  “sheath  boundary”,  the  problem  boundary,  or  throughout  all 
space.  Alternatively,  particles  can  be  initialized  with  a  user-generated  file.  Particles  are  tracked  for  a  specified 
amount  of  time,  with  the  timestep  automatically  subdivided  at  each  step  of  each  particle  to  maintain  accuracy.  The 
current  to  each  surface  element  of  the  spacecraft  is  recorded  for  further  processing.  The  charge  or  current  density 
created  by  the  particles  can  be  saved  for  use  in  solving  for  volume  potentials. 

III.  Plumes  in  Nascap-2k 

A  recent  enhancement  made  to  Nascap-2k  is  enhancement  of  the  ability  to  specify  ion  densities  due  to  an  ion 
source  (ion  thruster  or  plasma  contactor,  for  example)  and  use  them  in  the  computation  of  potentials  in  space. 
Potentials  due  to  plasma  sources  and  object  surfaces  are  needed  to  compute  contaminant  trajectories.  This  capability 
can  be  used  to  study  the  influence  of  spacecraft  surfaces  on  engine  plumes  and  interactions  between  plumes. 

Nascap-2k\  plume  capabilities  are  intended  to  complement  those  of  the  EPIC  (Electric  Propulsion  Interactions 
Code)  computer  code12  which,  like  Nascap-2k ,  was  developed  by  SAIC  for  NASA’s  SEE  program,  EPIC  is  a  fast- 
rnnning  code  that  models  plumes,  places  plumes  on  spacecraft,  displays  plume  densities  in  space,  and  calculates 
plume  fluxes  to  spacecraft  surfaces  and  resultant  surface  effects  (e.g.,  sputtering).  It  is  easy  to  modify  parameters  in 
EPIC ,  as  well  as  to  account  for  orbital  configuration  changes  (e,g.,  rotation  of  solar  arrays).  However,  EPIC  docs 
not  calculate  potentials  in  space,  and  thus  cannot  account  for  the  effects  of  surface  potentials  and  surface  sheaths, 
interactions  between  plumes,  or  self-consistent  flow  of  charge  exchange  ions  around  obstacles,  all  of  which  are 
intended  subjects  for  study  with  Nascctp-2L 

A  two-dimensional  (R 2)  map  of  the  ion  density  due  to  a  plasma  source  can  be  imported  directly  into  Nmeap-Ik. 
The  locations  and  directions  of  all  instances  of  the  plume  are  specified.  Nascap-2k  then  solves  a  nonlinear  Poisson 
equation  for  space  potentials  using  the  surface  potentials  as  boundary  conditions.  The  charge  density  consists  of  the 
known  ion  density  from  the  plume  specification  plus  a  barometric  electron  density: 
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where  0  and  n  are  the  plume  electron  temperature  and  reference  density  respectively,  p  is  the  total  charge  density  to 
be  used  in  Poisson’s  equation,  and  piDn  is  the  known  ion  density.  In  addition  to  the  specification  of  the  ion  density 
from  an  external  plume  map,  charge  exchange  ions,  such  as  those  generated  in  a  ion  thruster  plume,  can  be 
generated  and  tracked  within  Nascap-lk  to  achieve  a  self-consistent  solution13. 

This  capability  was  used  to  analyze  plume  interactions,  charge  exchange  ion  density  and  flow,  and  surface 
impingement  in  NEXT  multiple  thruster  tests14  conducted  at  NASA  Glenn  Research  Center  during  December  2005. 
The  test  configuration  includes  three  active  NEXT  thrusters  plus  a  dummy  thruster  in  the  fourth  position,  as  well  as 
a  full  complement  of  diagnostic  probes.  The  objective  of  the  test  was  to  identify  and  quantify  the  beneficial  and 
deleterious  Interactions  between  the  engines,  their  plumes,  and  their  controlling  electronics. 

Figure  2  shows  the  potentials  of  the  superposed  plumes  of  the  three  active  thrusters  on  a  plane  through  the  center 
of  two  of  the  thrusters.  This  is  used  as  an  initial  condition  for  the  calculations,  and  also  serves  as  a  point  of 
comparison  to  note  thruster  interaction  effects.  Except  immediately  adjacent  to  biased  surfaces,  the  potentials  are 
barometric,  /.e,  approximately  logarithmic  with  the  ion  density. 

Figure  3  shows  the  potentials  on  the  same  plane  as  in  Figure  2,  after  calculating  the  self-consistent  potentials  by 
generating  and  tracking  charge  exchange  ions  in  the  potentials  and  recomputing  potentials  until  a  self-consistent 
solution  is  reached.  Differences  of  note  between  the  self-consistent  potentials  (Figure  3)  and  the  superposed 
potentials  (Figure  2)  include  the  following: 

1.  Treatment  of  the  blockage  of  ions  by  the  thrusters  and  backing  plate.  The  blockage  of  ions  by  the  thrusters  and 
backing  plate  is  not  included  when  the  plumes  were  initially  generated  and  superposed,  resulting  in  substantial 
charge  exchange  ion  density  behind  the  blocking  plate.  The  self-consistent  calculation  shows  low  ion  density 
behind  the  plate  and  a  sensible  potential  and  density  structure  at  the  edge  of  the  plate. 

2.  A  charge  exchange  ion  stagnation  region  occurs  about  one  thruster  diameter  above  the  midpoint  between  the  two 
thrusters.  In  the  superposed  plume  structure,  charge  exchange  ions  are  accelerated  laterally  out  of  one  plume  and 
pass  unphysical ly  through  the  neighboring  plume,  whereas  in  the  self-consistent  structure  ions  are  blocked  by  the 
positive  potentials  of  the  neighboring  plume. 

3.  The  downstream  potentials  and  densities  are  higher  in  the  self-consistent  plume  than  in  the  superposed  plume, 
probably  because  of  inhibition  of  lateral  motion  by  neighboring  plumes. 
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Figure  2 *  Potentials  in  the  plumes  of  the  thruster  array  with  the  total  density  a  superposition  of  the  individual 
plumes,  used  as  initial  condition  for  the  calculation.  (Half  of  chamber  length  is  shown.)  Contours  are  in  a 
plane  through  the  center  of  two  active  thrusters*  Potential  and  density'  is  biased  toward  the  right  due  to  the 
dummy  thruster  at  the  left  rear. 
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Figure  3.  Self-consistent  potentials  in  the  plumes  of  the  thruster  array,  (Half  of  chamber  length  is  shown,) 
Contours  are  in  a  plane  through  the  center  of  two  active  thrusters.  Potential  and  density  is  biased  toward  the 
right  due  to  the  dummy  thruster  at  the  left  rear. 

IV.  Particle  in  Cell  Calculations 

Nascap-2k  has  primarily  been  used  to  model  quasi -static  phenomena.  However,  there  are  a  large  number  of 
physical  processes  of  interest  whose  timescales  require  a  dynamic  approach,  such  as  a  particle- in -cel  l  (PIC) 
technique.  Examples  of  such  processes  arc  breakdown  phenomena,  plasma  kinetics,  and  sheath  structure  about 
surfaces  with  potentials  that  change  on  a  timescale  comparable  to  the  lime  it  takes  an  ion  to  cross  the  sheath.  PIC 
techniques  can  also  be  used  to  address  problems  in  which  analytic  representations  of  the  environmental  currents  are 
inadequate,  such  as  in  a  spacecraft  wake  or  in  a  cavity.  A  PIC  technique,  in  which  the  ion  space  charge  density  is 
computed  from  macroparticles  tracked  from  the  boundary  of  the  computational  space  until  they  are  collected  or  exit 
the  computation  space,  was  successfully  used  to  model  the  CHAWS  experiment.  In  addition,  PIC  techniques  can  be 
useful  when  developing  analytic  models.  In  order  to  facilitate  these  modeling  techniques,  the  ability  to  perform 
various  types  of  PIC  calculations  was  built  into  Nascap-2k.  Recently,  additional  numeric  techniques  have  been 
added  to  make  two  types  of  PIC  calculations  more  flexible,  robust,  and  faster  Hybrid  PIC  and  full  PIC,  In  a  Hybrid 
PIC  calculation,  the  problem  is  initialized  by  creating  ion  macroparticles  throughout  the  grid  to  represent  a  constant 
particle  density.  The  ion  macroparticles  are  tracked  for  one  timestep  and  then  volume  potentials  are  computed  using 
the  resulting  ion  density  and  a  barometric  electron  density.  The  process  is  repeated  for  the  time  period  of  interest.  In 
a  full  PIC  calculation,  both  electron  and  ion  macroparticles  are  tracked  and  volume  potentials  are  computed  using 
the  resulting  plasma  density.  The  user  interface  was  also  revised  to  make  such  calculations  simpler  to  specify  and  to 
analyze. 
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Of  primary  interest  is  VLF  (about  1  to  20  kHz)  antennas  in  the  upper  ionosphere.  Plasma  waves  from  VLF 
antennas  with  such  frequencies  are  thought  to  interact  with  MeV  radiation  belt  electrons.  Such  an  antenna  would  be 
a  rod  several  inches  in  diameter  and  many  meters  long  and,  due  to  the  ease  of  electron  collection  by  positive  objects, 
would  be  nearly  always  at  negative  potential  relative  to  the  ambient  ionosphere.  Because  the  frequencies  of  interest 
are  comparable  to  the  ion  plasma  frequency,  the  sheath  structure  will  be  at  some  intermediate  state  between  the  “ion 
matrix”  or  “frozen  ion”  limit  (which  assumes  the  ions  are  stationary  and  contribute  ambient  ion  density  to  the  space 
charge)  and  the  equilibrium  space  charge  limit  (in  which  the  ions  assume  a  steady-state  space  charge  limited 
distribution  of  charge  and  current).  Thus,  calculation  of  the  sheath  structure  and  of  the  ion  collection  by  the  antenna 
requires  dynamic  (PIC)  treatment,  at  least  for  the  ions.  Nascap-2k  can  be  used  to  perform  all  four  simulations  of 
interest:  (1)  equilibrium  space  charge  sheath;  (2)  “frozen  ion”  sheath;  (3)  dynamic  PIC  ions  with  fluid  (Boltzmann  or 
barometric)  electrons  (Hybrid  PIC);  (4)  dynamic  PIC  ions  and  electrons  (Full  PIC). 

Two  of  these  recent  enhancements  to  Nascap-2k'$  PIC  capabilities  are  the  injection  of  macroparticles  from  the 
boundary  during  a  calculation  and  the  splitting  of  the  macroparticles. 

In  order  to  replace  mac roparl ides  that  are  collected  by  the  probe  or  escape  from  the  grid,  it  is  necessary  to 
periodically  inject  macroparticles  from  the  boundary.  This  allows  for  the  calculation  of  current  for  longer  time 
periods.  In  hybrid  PIC  calculations  without  boundary  injection,  the  low  field  region  near  the  boundary  of  the 
problem  develops  a  significant  negative  potential  due  to  ion  depletion.  Boundary  injection  keeps  these  potentials 
near  zero  by  replenishing  the  ions  that  have  been  collected  or  escaped. 

Particle  injection  is  implemented  with  an  injection  point  at  each  quarter-boundary-surface-elemcnt.  The  injected 
particle  has  charge  equal  to  the  plasma  thermal  current  times  the  area  times  the  time  interval,  and  velocity  equal 


2eT 

to  J- — —  ,  so  that  it  represents  the  inbound  half  of  the  plasma  density.  When  the  spacecraft  is  moving  through  the 
i  nm 


plasma  this  algorithm  is  modified  to  account  for  the  motion.  The  current  and  velocity  arc  computed  in  such  a  way 
that  the  density  contribution  of  the  injected  part  ides  varies  from  half  the  ion  density  for  a  stationary  object  to  the  full 
ion  density  for  a  high  mach  number  object. 

Closely  connected  with  boundary  injection  is  macroparticle  splitting.  There  are  two  major  reasons  for  splitting 
macroparticles:  one  physical  and  one  numeric.  Even  at  moderate  potentials,  thermal  effects  can  reduce  collected 
currents.  Some  particles  near  the  sheath  edge  have  enough  thermal  velocity  perpendicular  to  the  electric  field  that 
angular  momentum  conservation  prohibits  collection.  Particle  splitting  allows  for  a  representation  of  the  thermal 
distribution  in  the  initial  particle  distribution  and  in  particles  injected  from  the  boundary.  From  a  numeric  point  of 
view,  particle  splitting  can  be  used  to  keep  the  particle  weight  appropriate  to  the  grid  size  and  to  help  maintain  the 
smoothness  of  the  distribution.  A  large  particle  that  originated  in  an  outer  grid  is  split  so  that  it  becomes  distributed 
over  several  volume  elements  in  an  inner  grid. 

Macroparticles  may  be  split  when  they  enter  a  more  finely  resolved  region  or  when  they  arc  created  either  at  the 
boundary  or  throughout  the  volume  at  problem  initialization.  The  algorithm  used  is  as  follows; 

L  Particles  are  split  in  velocity  space  only.  Because  high-field  regions  are  often  of  interest,  spatial  splitting  would 
raise  problems  with  energy  conservation. 

2.  Each  particle  carries  a  temperature,  which  is  treated  as  isotropic.  The  fission  products  carry  half  the  temperature 
of  the  original  particle,  while  the  remaining  thermal  energy  appears  as  kinetic  energy  of  the  split  particles.  That 
each  macroparttcle  has  a  temperature  means  that  they  can  be  split  repeatedly  in  mid- flight,  each  time  the  particle 
enters  a  more  Finely  resolved  region. 

3.  For  splitting  purposes  the  Z-axis  is  defined  to  be  along  the  direction  of  the  particle  velocity,  the  X-axis  randomly 
chosen  in  the  plane  normal  to  Z,  and  the  Y-axis  mutually  perpendicular, 

4.  Particles  are  split  into  two  or  three  particles  along  each  axis,  except  that  a  particle  is  not  split  along  the  Z- 
direction  if  the  kinetic  energy  exceeds  the  thermal  energy.  Not  splitting  along  Z  helps  avoid  particle 
proliferation,  but  makes  an  error  by  not  preserving  the  original  particle  temperature  along  Z.  Eight,  nine,  or 
twenty-seven  new  particles  result. 

5.  Particle  velocity  is  assumed  to  be  acquired  by  acceleration  rather  than  actual  drift  (i.e.,  spacecraft  velocity).  If 
there  is  actual  drift,  then  the  drift  velocity  is  removed  before  splitting  the  particle  and  added  back  after, 

6.  If  the  particle  is  split  by  two  along  the  X  or  Y  axis,  the  new  velocity  is  ±0,707  ,/r/m  .  Along  the  Z  axis  the 


muG  J 


y 

velocity  increment  is  calculated  as  if  the  temperature  were  T  -  2muQ 
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7.  If  the  particle  is  split  by  three  along  the  X  or  Y  axis,  there  is  a  zero-velocity  central  panicle  and  two  “probe” 
particles  with  velocity  is  ±0.866^/7/01  .  Along  the  Z  axis  the  velocity  increment  is  calculated  as  if  the 

/  \ 

2 

temperature  were  T-2muQ 

8.  The  new  particles  have  the  same  properties  as  the  original  particle  except  for  velocity,  weight  (charge),  and 
temperature. 

We  use  Nascap-2k  to  compute  the  current  collected  by  a  10  cm  radius,  -100  V  sphere  in  a  I0m  m  \  1  eV,  H+ 
plasma.  The  sphere,  represented  by  the  32  surface  quasisphere  shown  in  Figure  4,  is  embedded  in  a  2.4  m  cubic  grid. 
The  plasma  density  is  computed  using  the  Hybrid  PIC  charge  density  model  with  timesteps  of  5  x  10 7  s.  The 
collection  of  ions  by  the  probe  and  the  loss  of  ions  out  the  sides  are  monitored,  and  the  final  potential  and  particle 
configurations  are  inspected. 


Figure  4.  Sphere  object  used  in  example. 

The  sphere  is  a  useful  test  case  not  only  because  the  equilibrium  current  can  be  analytically  computed,  but 
because  the  symmetry  provides  a  good  test  of  the  ability  of  the  splitting  to  properly  represent  the  angular  momentum 
barrier  An  analytic  formulation  gives  the  equilibrium  sheath  radius  of  73  cm  for  a  collected  current  of  41 
microamperes.  Presheath  enhancement1 5  gives  another  factor  of  1.45  for  59  microamperes.  Presheath  enhancement 
is  the  focusing  of  current  outside  the  sheath.  The  PIC  approach  automatically  includes  the  presheath  enhancement. 
The  current  leaving  the  grid  is  the  surface  area  of  the  grid  times  the  one-sided  plasma  thermal  current  for  a  total  of 
216  microamperes. 

Figure  5  compares  the  collected,  leaving,  and  injected  current  for  a  calculation  with  particle  splitting  and 
boundary  injection  and  a  calculation  with  neither.  The  collected  current  rapidly  rises  to  a  sustained  value  of  about  60 
microamperes.  Without  splitting  the  current  fluctuates  widely.  With  no  splitting,  there  is  no  current  leaving  the  grid 
as  the  electric  field  points  toward  the  grid  center  and  there  is  no  thermal  motion.  With  splitting  and  injection,  the 
escaping  and  injected  currents  arc  much  larger  than  the  collected  current,  averaging  160  to  170  microamperes.  The 
injected  current  is,  on  average,  slightly  greater  than  the  escaping  current,  and  far  less  than  the  sum  of  the  collected 
and  escaping  currents. 

Figure  6  and  Figure  7  show  the  potential  contours  for  the  two  calculations  after  50  microseconds.  The  potentials 
are  spherical.  While  the  potentials  with  no  splitting  or  boundary  injection  look  smoother,  after  50  microseconds  they 
are  extending  further  from  the  sphere  as  the  total  ion  charge  in  the  grid  is  depleted.  The  potential  fluctuations  in  the 
case  with  splitting  and  boundary  injection  are  all  well  under  the  plasma  temperature  of  1  cV. 

For  comparison,  Figure  8  and  Figure  5  show  results  for  the  case  with  splitting  but  no  particle  injection.  In  this 
case,  between  the  sheath  edge  and  the  grid  boundary,  the  potential  is  more  negative  than  -0.3  volts.  This  is  because 
the  ion  population  is  depleted  due  to  both  collection  of  ions  by  the  probe  and  by  escape  of  ions  through  the  grid 
boundary  due  to  their  thermal  motion.  Significant  potential  drop  can  be  caused  by  a  relatively  small  depletion,  so 
that  at  50  microseconds  the  collected  current  is  not  significantly  reduced. 
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These  improvements  come  at  a  considerable  cost  in  computing  time  and  disk  space.  From  an  initial  particle 
count  of  0.3  million  for  the  case  with  no  splitting  and  no  injection,  the  particle  count  decreases  to  0J  million  at  50 
microseconds;  while  the  run  with  particle  splitting  and  boundary  injection  starts  at  2.4  million  particles  and  increases 
to  8.4  million  at  50  microseconds. 

Figure  9  shows  a  quadrant  of  particles  after  nine  microseconds.  The  particles  are  initially  unsplit,  and  the  ions 
arc  split  into  nine  or  twenty-seven  particles  on  entering  a  finer  grid.  In  Figure  9  ions  can  be  seen  entering  Grid  2 
from  Grid  I  at  the  top  and  right.  Because  these  ions  are  moving  slowly,  they  are  split  both  along  and  normal  to  their 
motion  direction.  By  this  time  all  ions  originating  in  Grid  3  have  been  “eaten"  by  the  sphere,  so  that  the  cloud  of 
ions  currently  in  Grid  3  have  entered  from  Grid  2  and  been  split.  Those  that  entered  most  recently  were  already 
drifting  significantly  when  they  were  split,  and  were  thus  split  only  normal  to  their  direction  of  motion. 
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Figure  5.  Collected  current  (left  scale),  escaping  current  (right  scale),  and  injected  current  (right  scale) 
running  problem  with  boundary  injection. 
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Figure  6.  Potential  contours  after  50  microseconds  for  Hybrid  PIC  calculation  with  no  particle  splitting  and 
no  boundary  injection. 


Figure  7.  Potential  contours  after  50  microseconds  for  Hybrid  PIC  calculation  with  particle  splitting  and 
boundary  injection. 
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Figure  8.  Potential  contours  after  50  microseconds  for  Hybrid  PIC  calculation  with  particles  split  upon 
creation  and  upon  entering  grid  with  finer  resolution,  but  no  boundary  injection. 


Figure  9.  Particle  positions  after  nine  microseconds  when  particles  are  split  on  entering  a  more  finely  resolved 
grid. 
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The  thermal  distribution  of  the  ions  plays  an  important  role  in  the  formation  of  the  wake  structure.  To  illustrate 
this,  we  examine  the  wake  structure  with  and  without  splitting.  Figure  10  shows  a  hybrid  PIC  calculation  (with  no 
particle  splitting)  of  O  flow  past  an  uncharged  sphere  moving  at  LEO  velocity  in  the  (1,1,0)  direction.  Clearly  seen 
are  the  negative  potentials  tn  the  object  wake  and  the  boundary  between  the  injected  particles  (diagonal  pattern)  and 
the  original  particles  (square  pattern).  Potential  fluctuations  on  the  order  of  0.05  volts  are  seen  in  the  first  subdivided 
grid  where  it  is  populated  by  outer  grid  ions.  Figure  1 1  shows  a  similar  calculation,  now  with  the  sphere  at  -100  V. 
After  1 36  microseconds  ions  are  focused  in  the  wake  sufficiently  strongly  to  create  positive  potentials  near  the  ram 
energy.  (Note  that  this  rather  striking  result  is  probably  associated  with  the  object  size  being  comparable  to  the 
plasma  debye  length.) 

Figure  12  shows  the  charged  sphere  calculation,  now  with  particles  split  on  entering  a  refined  grid.  While  the 
general  character  of  the  result  is  the  same,  the  potentials  are  much  smoother  and  now  show  compression  of  the 
sheath  on  the  ram  side.  Figure  13  shows  the  current  to  the  sphere.  After  an  early  peak  to  nearly  16  microamperes, 
the  current  settles  down  to  a  value  of  less  than  7  microamperes,  comfortably  less  than  the  orbit-limited  value  of  8 
microamperes.  Of  course,  this  improved  Fidelity  comes  at  a  price,  with  a  final  particle  count  in  excess  of  two  million 
in  Figure  12,  versus,  under  0,4  million  in  Figure  ft. 


Figure  10.  Potentials  and  ion  (0+)  macroparticles  after  80  microseconds  for  an  uncharged  sphere  moving  in 
the  (1,1,0)  direction.  (No  splitting  of  macroparticles.) 
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Figure  11,  Potentials  and  ion  (O*)  macroparticles  after  136  microseconds  for  -100  V  sphere  moving  in  the 
(1,1,0)  direction.  (No  splitting  of  macroparticles.) 


Figure  12*  Potentials  and  ion  (O*)  macroparticles  after  160  microseconds  for  a  sphere  charged  to  -100  V 
moving  in  the  (1^1,0)  direction.  Particles  split  on  entering  refined  grid. 
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Another  new  capability  is  orbit  averaging.  Previously  Nascap-2k  could  perform  two  kinds  of  particle  in  cell 
calculations.  In  one  kind,  used  to  model  the  CHAWS  experiment,  which  involved  current  collection  by  a  high 
potential  in  the  wake  region  of  a  disk  shaped  spacecraft  (WSF),  macroparticles  arc  created  at  the  boundary  of  the 
computational  space  and  arc  tracked  until  they  either  are  collected  by  the  object  or  leave  the  computational  space 
These  macroparticles  carry  current  and  share  the  current  times  the  sub-step  time  in  each  volume  clement  each  sub¬ 
step.  After  all  the  macroparticles  are  tracked,  potentials  are  computed  using  the  ion  charge  deposited  by  the  tracked 
macroparticles.  The  full  trajectory  tracking  is  iterated  with  the  potential  solution  until  a  self-consistent  solution  is 
found.  Often  mixing  of  the  present  and  previous  solutions  is  necessary  to  reach  self-consistency.  The  other  type  of 
PIC  calculation  is  the  traditional  approach  discussed  above,  in  which  the  macroparticles  are  only  tracked  fora  single 
timestep  and  their  entire  charge  is  deposited  on  the  computational  nodes  at  the  end  of  the  timestep.  This  approach  is 
needed  for  time-dependent  phenomena. 

Orbit  averaging  is  a  variation  of  the  standard  PIC  approach  in  which  a  portion  of  the  macroparticie's  charge  is 
deposited  to  the  gnd  each  substep.  This  approach  is  appropriate  when  the  timescales  of  the  phenomena  of  interest 
are  long  compared  with  the  time  it  takes  a  panicle  to  transit  a  volume  element.  The  timestep  can  be  appropriate  to 
the  phenomena  of  interest,  while  the  substeps  are  such  that  the  particles  travel  only  a  fraction  of  mesh  unit.  This 
approach  can  give  a  much  smoother  solution  as  seen  in  Figure  14,  A  variation  will  be  appropriate  for  the  tracking  of 
electrons  when  the  timestep  is  set  by  the  ion  transit  time. 
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Figure  14.  Current  to  sphere  with  standard  PIC  approach  and  with  orbit  averaging. 

The  final  improvement  we  discuss  is  the  ability  to  compute  the  change  in  surface  potentials,  charging,  due  to 
tracked  ions  from  a  Hybrid  PIC  or  full  PIC  calculation.  This  capability  makes  it  possible  to  calculate  spacecraft 
floating  potential  and  differential  charging  for  cases  where  analytic  approximations  to  plasma  currents  are  not 
adequate,  for  example  with  a  dynamic  bias  or  a  concave  geometry. 

When  this  technique  is  usedt  ions  arc  tracked  and  the  net  surface  current  is  the  sum  of  the  tracked  ion  current  and 
an  analytic  electron  current.  The  electron  current  is  a  function  of  the  surface  area.  A,  the  electron  thermal  current,  jlh, 
the  surface  potential,  +,  and  the  plasma  temperature,  0, 

Aith  «p(+/9) 

i^H) 

Wc  compare  the  analytic  and  numeric  calculations  of  the  discharge  of  a  0, 1-m  radius,  -100  V  sphere  in  a  1010  m 
\  1  eV,  hydrogen  plasma.  The  Nascap-2k  calculation  is  done  using  a  32  surface  quasisphere  shown  in  Figure  4.  The 
plasma  density  is  computed  using  the  Hybrid  PIC  charge  density  model  with  timesteps  of  5  x  10 7  s. 

The  change  in  potential  of  a  sphere  of  radius  r  in  vacuum  due  to  an  incident  current  I  during  a  timestep  of  length 

t  is  given  by  A+  =  ~  where  C  =  47t£0r.  Figure  15  compares  the  potential  computed  by  Nascap-2k  with  the 
C 

potential  computed  from  the  capacitance  of  a  sphere  of  radius  0 A  m  and  the  tracked  current.  The  figure  also  shows 
the  potential  computed  for  the  capacitance  of  a  sphere  of  radius  0,09285  m. 

With  the  analytic  electron  current,  the  total  current  drops  to  zero  as  the  potential  becomes  slightly  positive,  and 
the  potential  is  held  near  zero  by  the  incident  electrons.  The  ion  current  remains  at  a  fairly  constant  rate  due  to  the 
inertia  of  the  ions.  Since  the  ion  current  is,  on  average,  slightly  less  than  the  raw  electron  thermal  current,  the 
potential,  on  average,  is  slightly  negative. 
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Figure  15.  Current  to  sphere  and  resulting  potentials  for  sel  ^consistently  computed  discharge. 


V,  Conclusion 

The  recent  enhancements  of  Nascap-2k  will  allow  users  to  analyze  a  wider  variety  of  spacecraft- pi  asms 
interactions.  The  improvements  to  the  ability  to  specify  ion  densities  due  to  an  ion  source  and  use  them  in  the 
computation  of  potentials  in  space  allows  for  modeling  of  thrusters  plasmas.  The  new  PIC  capabilities  will  allow 
users  to  address  fully  dynamic  problems  in  addition  to  quasistatic  interactions  for  which  it  is  already  widely  used* 
The  capabilities  described  here  will  be  available  in  the  next  code  release  later  this  year.  Na$cap-2k  is  distributed 
through  the  SEE  program  web  site,  http://see.  msfe.  nasa.gov. 
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